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Impulses efficiently propagate into nominally dry granular beds and backscatter from buried inclu- 
sions in such beds may be potentially exploited to image shallow buried objects (SBOs). However, 
reliable imaging of SBOs requires "cleaning up" of surface vibrations, and, in addition to 3D particle 
dynamics simulations, a phenomenological model to parameterize the bed surface may be useful for 
held applications. We introduce a ID mean-field-like toy model with two parameters that allows one 
to model surface vibrations, is consistent with experiments in a granular bed, and can help estimate 
the approximate signal transmission properties of the bed. 
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Introduction: The imaging of shallow buried objects in 
a complex medium, e.g., nominally dry soil, is a difficult 
problem that has seen limited progress Such imaging 
is of relevance in connection with locating antipersonnel 
land mines, in archaeology, land surveying and in other 
applications. It has been shown that gentle mechani- 
cal impulses Q can be used to detect buried objects at 
depths of a meter or so in nominally dry sand beds. De- 
tailed 3D simulations establish that nonlinear pulse prop- 
agation in 3D beds is a quasi- ID process 0]; normally 
incident pulses travel as weakly dispersive energy bun- 
dles and become more and more ID-like with increase in 
area over which the impulse is generated. It would be of 
interest to rapidly generate images of buried backscatter- 
ers by exploiting the information contained in the time- 
dependent surface vibrations in granular beds 0- 

To accomplish such imaging, it is necessary to probe 
some global parameter that contains coarse-grained in- 
formation about grain dynamics at the bed surface. We 
study the space averaged, time evolution of the time in- 
tegrated kinetic energies of the surface grains in ideal- 
ized sand beds. Several groups have measured impulse 
backscattered signals at the surface in empty beds and in 
beds with some buried object Q. Newtonian dynamics 
based 3D simulations of impulse backscattering in ideal- 
ized beds have been carried out by Sen et al. |3J. These 
authors used the velocity- Verlet algorithm Q to integrate 
the equations of motion, and their data are consistent 
with the available experiments. Thus, a consensus on the 
spatio-temporal behavior of impulse backscattered data 
is beginning to emerge. 

It is apparent from the existing work that the energy 
imparted by the impulse penetrates into the system. The 
spread of the energy in a given x-y plane at a given depth 
depends on the packing in the system. There is impulse 
backscattering at every granular contact. If one measures 
the amount of backscattered energy at the bed surface, 
the energy density at the bed surface rapidly depletes 



TABLE I: Description of the energy transfer process in the 
ID model. The arrows indicate the direction in which the 
layers (lyrl etc.) will interact in the next time (t) step. The 
exchange case has been depicted for convenience. e hs (t) is the 
backscattered energy received back at the surface. 
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after the initiation of the impulse, and then rises as a 
function of time. The backscattering from the shallow 
layers is significant. The amount of backscattering from 
the deeper layers does get weaker, and eventually dies 
out. The dissipative properties of the bed play an im- 
portant role in the attenuation of the impulse. 

We contend that, at least for the purposes of held ap- 
plications, it is desirable to explore a tractable ID toy 
model that can capture the critical results of the 3D sim- 
ulations. In this Communication, we propose a two pa- 
rameter model to describe impulse backscattering at the 
surface of a granular assembly. It may be necessary to 
use more parameters (e.g., to include information about 
the area across which the impulse is imparted), if one 
is looking for detailed agreement with experiments. The 
physics is similar in spirit to that of mean-held theories. 

Phenomenological Model: We define a vertical align- 
ment of layers, where each layer can be thought of as a 
mass 0. At time t = 0, we set initial energy E = 1 for 
layer one and zero for the rest. At t — 1, the first layer in 
the vertical chain transfers p (< 1) of the impulse energy 
to the second layer, and retains (1 — p). At subsequent 
times, the impulse will propagate in the same fashion, at 



every step, all the way down the chain (see Tabled . Each 
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push the preceding layer in the opposite direction in the 
following time step. Since the phase reverses every time 
step, they will interact, alternately, with layers above and 
below, in alternate time steps, thereby introducing sig- 
nificant backscattering into the problem. Note also that 
successive layers will be in opposite phases at any time, 
assuring interaction between them only in alternate time 
steps (Table HJ|. 

We model the interaction between two adjacent layers 
in our simulations in the following two ways: (i) equipar- 
tition case, and (ii) exchange case. In the equipartition 
case, the two interacting layers will come away from the 
interaction with equal amounts of energy; we add up the 
individual energies of the two layers and divide the sum 
equally between them. In the exchange case, we let the 
layers exchange their energies; the two interacting lay- 
ers, after the interaction, come away with the energy of 
the other. The equipartition case can be viewed as one 
that leads to ergodic-like behavior, where we assume that 
the two adjacent layers get compressed to the same ex- 
tent during the interaction, and the potential energy of 
compression gets converted back to the respective kinetic 
energies, which will now be half of the total energy of the 
two layers. The equipartition ansatz negates the symme- 
try breaking introduced by p ^ 0.5. The exchange model 
captures the essence of nonlinear impulse propagation in 
which an impulse travels as a perfect solitary wave in a 
ID chain of elastic grains 0, and as a weakly dispersive 
energy bundle in 3D beds we model the situation 
where two energy bundles, traveling in opposite direc- 
tions, go through each other without distortion. In real 
systems, one would expect that both the equipartition 
and exchange behaviors would be present, and such an 
extension of our study will be reported elsewhere |2j. 

We monitor the energy transfer at the surface in our 
model analysis. The first layer, in its negative phase (we 
assume positive phases to point down the chain), will 
transmit p fraction of its energy to the surface and retain 
(1 — p) fraction to itself; the surface does not transfer 
any energy back to the first layer. At the bottom of the 
chain, we let the last layer lose p fraction of its energy 
in its positive phase and retain (1 — p); the lost energy 
is presumed to travel further down in a similar fashion. 
These boundary conditions do not, in any way, affect our 
final results. We have verified our results with longer 
chains and there are no qualitative changes (see further 
discussion below); we have, therefore, employed a 40 layer 
long system for our studies. 

Some Analytical results: For the exchange model, the 
sequence of backscattered energy packets that arrive at 
the surface, e bs (t), can be worked out (c/ Table H] see 
Fig. da)): 

p(l - P ),p 2 (l -p),p 3 (l -p), ■ ■ ■ ,p(l-pf,p 2 (l~pf, 
p 3 (l-pf, ■ ■ ■ ,p(l-p)\p 2 (l - p)\ P 3 (l - pf, ■ ■ ■ 

and so on. Since p < 1, the sum of the first subsequence 
(each subsequence is separated by • • •) yields p, and, sim- 
ilarly, the sum of the second subsequence yields p(l — p) 2 



and so on. Thus, the entire sequence of energy packets 
received back at the surface can be summed to obtain, 

with E bs (t) respresenting J2t'=o e&s (*')> 

E bs (t = oo) = p[l + (1 - pf + (1 - pf + • • ■] - 1/(2 -p) 

We see that E bs (t = oo) will always be > 1/2; also, 
E bs (t — oo) is higher if p is higher. If p < 0.5, the suc- 
cessive packets of e bs (t) in each subsequence is going to 
decrease fast, and, when the next subsequence involving 
a higher power of (1 — p) arrives, a jump in E bs (t) takes 
place. In the other case, when p > 0.5, there is going 
to be a faster increase of E bs (t) within each subsequence 
itself, than brought in by the arrival of the next subse- 
quence. The above analysis shows that a plot of E bs {t) 
against time is going to be characterized by a sequence of 
steps and plateaus, as is seen from Fig. da), first inset. 
The subsequences can be understood as due to the soli- 
tary wave like propagation 0. where each packet trav- 
els down the chain undiminished till it loses energy (by 
p fraction) at the two ends. Note that the sequences are 
geometrically decreasing, and, as may be gathered from 
the plot (top smaller inset), characterized by an exponen- 
tial decrease in packet sizes. As a consequence, E bs (t) 
is characterized approximately by logarithmic growth in 
the subsequences (bottom smaller inset), with the added 
feature of a saturation plateau seen in the earlier inset. 

An analysis like the above is not possible for the 
equipartition case because of the complicated algebra re- 
sulting from distributing half of the total energy of the 
two particles to each, after the interaction. Neverthe- 
less, numerical results indicate that a similar result, of 
E bs (t — oo) being always > 1/2, holds in this case as 
well. Also, there are no steps and plateaus in this case, 
except for the ones introduced by the fact that e bs (t) 
arrives in alternate time steps (Fig. ^b); see further dis- 
cussion below). 

Numerical Analysis: Soil is a highly heterogeneous 
medium, and impulse propagation and scattering prop- 
erties vary much spatially, decreasing and/or increas- 
ing with successive layers. Therefore, we investigate the 
above models for different scenarios characterized by dif- 
ferent distributions of p. A reasonable model that can be 
chosen, for qualitative comparison of results for different 
p's, is that of a logarithmic growth pattern in E bs (t). 
This is best illustrated for p = 0.25 (for all layers) in 
the equipartition case (Fig. ^b)); depending on the dis- 
tribution and values of p, the growth can be higher or 
lower than an average logarithmic growth in the initial 
stages. In Fig.^b), it is seen that e bs (t) = in alternate 
steps; also, the packets arrive in exponentially decreas- 
ing quanta in alternate steps. Such decay in e (t) is 
expected in view of the fact that most of the backscat- 
tered energy, at early times, is from the first few layers, 
where, by construction, the energy will be more. The 
deeper layers, the backscattered energy of which reaches 
the surface at later times, contribute to progressively 
higher order backscattering. The first inset to the fig- 
ure shows E bs (t), growing rapidly in time at early times, 
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FIG. 1: The two cases of exchange and equipartition axe shown in (a) and (b) respectively, where p is constant across the layers. 
The two smaller insets in (a) show the exponential decay of e bs (t) within a subsequence (top; y-axis logarithmic), and the 
associated logarithmic growth in E ba (t) (bottom; x-axis logarithmic); the solid lines are the corresponding curve fitting lines. 
Smaller insets in (b) show the logarithmic growth of E b "(t), in the initial stages, in the equipartition case (bottom inset uses the 
same data as the top inset, plotted with logarithmic a>axis). Available 3D simulation results || show favorable comparisons 
with these results. 



followed by progressive slowing down. The two smaller 
insets in this figure shows, as a zoom-in, the region close 
to the origin, of the earlier inset, to show the form of the 
(logarithmic) growth in E bs (t); the bottom inset plots 
the data with logarithmic abscissa for explicit confirma- 
tion. This pattern has been seen in 3D simulations of 
impulse propagation in sand beds @, and is consistent 
with experimental results 0. 

The maximum attained value of E bs (t), referred to as 
Emaxi 1S high compared to the 3D simulation results. 
Thus, the ID toy model, with its restrictions in the avail- 
able energy channels tends to backscatter much more en- 
ergy to the surface than is typical of 3D beds. There are 
two ways to improve the model to better mimic the prop- 
erties of 3D beds, (i) by incorporating restitution between 
layers, and (ii) by varying p appropriately as a function 
of position to approximately account for changes in layer 
compression as a function of depth, inhomogeneities in 
the medium etc. We discuss the effects of these decora- 
tions in our model below. 

Introducing a tunable parameter, q (< 1), to account 
for restitutional losses at each interaction (as also the 
spreading away of energy in 3D), which works to dissipate 
energy, allows us to tune E^ ax ; at each step, a fraction 
(= q) of the total energy of the two interacting masses 
is removed. It was observed (see Fig.EJa) and (b)) that 
E^ ax decreases exponentially with q; it drops by an order 
of magnitude as q is varied from 0.1 to 0.5. 

Now, we consider the richness of the model by explor- 



ing various distributions of p and their effect on the pat- 
terns of e bs (t) . We will only discuss the results using the 
plots of E bs (t) which arc, firstly, easier analysed, because 
of the logarithmic growth pattern in the initial stages. 
Secondly, E bs (t), being an integrated quantity, is, per- 
haps, more suited for comparison with experimental re- 
sults. Hence, we will plot these graphs on a semi-log basis 
for the study, with time axis being taken logarithmically. 

Distributions ofp and corresponding patterns in growth 
of E bs (t) : With a uniform distribution of p's (p constant 
across the layers), analysis of the equipartition model 
shows that the growth is, qualitatively, symmetric with 
respect to p = 0.5; E^ ax is the highest for p = 0.5, 
and falls off when p ^ 0.5. Clearly, the equipartition 
model is working best when the initial impulse itself is 
transmitted down the chain with equal sharing of the en- 
ergy between the leading layer and the next one in the 
chain. On the other hand, the curves are also symmet- 
ric, on either side of p = 0.5, around the mid-value of 
the respective ranges, i.e. around 0.25 and 0.75. The lat- 
ter symmetry is with respect to the initial growth rates. 
The curves are concave (lesser than average logarithmic 
growth rate), in the initial stages, on the semi-log plots 
for p < 0.25 (p > 0.75), and convex (larger than average 
logarithmic growth rate) for p > 0.25 (p < 0.75). This 
symmetry is natural because of the symmetry introduced 
by the exchange of roles between p and (1 — p) as p is var- 
ied. We get the best (approximate) logarithmic growth, 
in the initial stages, for p = 0.25 (p = 0.75). 
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FIG. 2: Growth patterns of E bs (t) with an inverse square law 
type variation in p. In (a) and (b), equipartition and exchange 
models, respectively, are shown for p decreasing towards the 
bottom; in (c) and (d), similarly, for p increasing towards the 
bottom. Note that the y-ranges in (c) and (d) are extremely 
small compared to those in (a) and (b). The effect of intro- 
ducing a dissipation parameter, q, is shown in (a) and (b). 

In the exchange case, there is no such symmetry. It is 
a monotonous behavior with E hs (t) requiring many more 
cycles to saturate, the lower the p value is, and saturating 
faster, the higher the p value is. In the latter case, if 
the p value is sufficiently high, E^ ax can be attained 
almost in the first cycle. This is understandable since 
the undispersed energy bundles that travel up and down 
the chain, but for the losses at the two ends, are smaller 
(larger) in size when p is smaller (larger) . 

We may let p decrease towards the bottom and vice 
versa. In the equipartition case, if we impose a linear de- 
crease (increase), in the initial stages, the graphs are con- 
vex (concave) and the variation in the amount of decrease 
(increase) shows variation in convexity (concavity). This 
is due to the longer presence of energy in the upper lay- 
ers, with higher (lower) p values . Similarly, for the ex- 
change model as well, steeper initial growth is achieved 



for higher p's in the top layers and vice versa. We may 
impose other rates of decrease (increase) like exponen- 
tial, inverse square law etc. but the basic charateristic 
of the convexity (concavity)/steepness of the equiparti- 
tion/ exchange graphs being dictated by the presence of 
larger (smaller) p's in the upper layers will remain the 
same; only it will be more pronounced with larger rates of 
decrease (increase) like exponential or inverse square law 
etc. In Figs.UJa) and (b), we have shown results for a case 
where p's decrease in inverse square law fashion towards 
the bottom for the equipartition and exchange models, 
respectively; note the high values of E^ ax . Figs. [2fc) 
and (d) show the results, similarly, for a case where the 
p's increase towards the bottom; note that, if the p's in 
the upper layers are extremely small (of the order of, say, 
10 -3 ), E^ ax is going to be extremely small (< 10~ 3 ) in 
both the models. 

If we let the p's vary randomly, we get a combination 
of the effects discussed earlier and the model is sensitive 
to different realizations of the random distribution of p's, 
a desirable feature for a simple phenomenological model 
for a highly heterogeneous medium like soil. We have 
also checked our results with longer chains. It is clear 
that the cyclicity in the exchange model will change with 
the length of the chain; the cyclicity is twice the length of 
the chain because we are only monitoring the energy at 
one end. This introduces, for example, a longer plateau 
length at each step of the staircase, in certain situations. 
However, E^ ax remains the same in the case of a con- 
stant p across the layers. In the equipartition case, E b f ax 
does change slightly; however, the convexity (concavity) 
patterns in the plots, which depend on the p values, are 
not changed, only they get more accentuated. 

In summary, it can be stated that our ID model is ro- 
bust, and, at the same time, sensitive to changes in the 
distribution of p across the layers. The model qualita- 
tively reproduces the time evolution of average energy of 
the surface grains of a complex, heterogeneous medium 
like granular beds, following the generation of a normal 
impulse at the bed surface. The model promises to be a 
useful tool for further research on this important topic. 
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